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Abstract 

By means of a simple model system, the total volume fluctuations of a tapped granular material 
in the steady state are studied. In the limit of a system with a large number of particles, they 
are found to be Gaussian distributed, and explicit expressions for the average and the variance are 
provided. Experimental and molecular dynamics results are analyzed and qualitatively compared 
with the model predictions. The relevance of considering open or closed systems is discussed, as 
well as the meaning and properties of the Edwards compactivity and the effective (configurational) 
temperature introduced by some authors. Finally, the linear response to a change in the vibration 
intensity is also investigated. A KWW decay of the volume response function is clearly identified. 
This seems to confirm some kind of similarity between externally excited granular systems and 
structural glasses. 

PACS numbers: 47.50. Cc, 05.50.+q, 81.05.Rm 



1 



I. INTRODUCTION 



Compaction can be roughly defined as the density relaxation of a loosely packed system 
of many grains under mechanical tapping or vibration l|, l2|. Both the time relaxation of the 
system and the steady state eventually reached in the long time limit, exhibit interesting 
properties, and have received a lot of attention during the last years. The description of 
this phenomenon goes well beyond the scope of equilibrium statistical physics and even of 
usual non-equilibrium theories. The time evolution of the system corresponds to a consec- 
utive series of mechanically stable states in which all the particles are at rest, at least at a 
macroscopic level of description. Therefore, the system is not thermal, in the sense that the 
fluctuations of the particles around their average positions play no role in characterizing the 
state of the system. In fact, in the description of the phenomenology of compaction, time is 
usually measured in number of taps that, of course, has nothing to do with the real motion 
of the particles. 

Edwards and collaborators formulated some years ago a theory providing a 
thermodynamic-like description of granular systems under certain circumstances. The main 
hypothesis of the theory is that for externally perturbed powders in a steady state, all the 
mechanically stable (metastable) configurations of a granular assembly occupying the same 
volume have the same probability. Therefore, it is assumed that the volume plays here a 
role analogous to that of the energy in usual thermal systems. The parameter conjugated 
to the volume, similar to the thermal temperature, was named 'compactivity'. 

In principle, the conditions under which the theory should apply are fulfilled in the 
compaction process. For this reason, the predictions of the Edwards theory have been 
compared with results obtained in experiments M, a ol, numerical simulations 3, M, and 
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simple models [9|, [lOl, lUl, 112|, 113|, [IJ, [151, [161, 117| of compaction. Although in general a fairly 



good agreement is claimed, at least in the limit of weak tapping, an analysis of the situation 
shows that there are many fundamental aspects that can not be considered as reasonably 
well established. A few of them will be addressed, although not solved with generality, here: 

1. Given that many of the experimental determinations of the compactivity are based in 
the measurements of density fluctuations, special care must be given to differentiate 
between fluctuations in the number of particles in a system of constant volume, fluc- 
tuations of volume in a system with a given number of particles, and fluctuations of 
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density in a system whose volume and number of particle change in time. All these 
fluctuations are not equivalent, in principle [l7 |. 

2. Some authors have reformulated the thermodynamic theory of powders using the (con- 
figurational) energy as the variable characterizing the macroscopic state of the system 
and an effective temperature, i.e. they use a formulation much closer to ordinary sta- 
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ations does not 



Il8 |. For instance. 



tistical physics. The equivalence or relationship between bot 
seem to be clear, although some ideas have been presented 
are the compactivity and the effective temperature independent variables or there is 
a relationship between them? 



3. The thermodynamic theory of powders was formulated for complete systems, and not 
as a local theory. In this context, it could be thought that a column of tapped granular 
matter in the steady state is characterized by a unique compactivity, aside from the 
boundary layers, even though the system is inhomogeneous. Another possibility is to 
consider that the system can be divided into horizontal layers, each of them having 
a different compactivity (or effective temperature). In fact, this local version of the 
theory has been often used in describing experiments and particle simulations. Which 
are the physical meaning and the implications of a gradient of the compactivity, or 
effective temperature, are also open questions. 

Another important aspect of compaction refers to the response to an external perturba- 
tion, e.g. a change of the vibration intensity. It happens that many of the peculiar dynamic 
behaviors exhibited by granular materials under tapping are similar to those shown by con- 
ventional structural glasses. For instance, both kind of systems present slow relaxation, 
annealing properties, and hysteresis effects. A necessary first step in order to understand 
these issues is to study the response of the system to a small perturbation of the conditions 
defining its steady state (linear response). Of course, the obvious external parameter to vary 
in the compaction experiment is the vibration intensity characterizing the tapping process, 
and the quantity whose response to analyze is the volume of the system. In the process, all 
the system, with a constant number of particles, is considered. 

An additional reason to study linear response in granular systems, is that fluctuation- 
dissipation relations obtained in molecular systems in the context of linear response theory, 
are being used to test the thermodynamic theory of granular matter {3, H]. Nevertheless, 



given the non thermal character of granular systems, even if the validity of the thermody- 
namic theory as proposed by Edwards is assumed, there is no reason a priori to expect the 
fluctuation-dissipation relations to apply. 

The aim here is to investigate the above mentioned issues by means of a very simple 
model for compaction [itI. The model dynamics is defined in such a way that it conserves 
the number of particles in the system, i.e. it is a closed system. Its volume changes due 
to the creation and destruction of empty sites or holes in the lattice over which the system 
is defined. The model is consistent with the Edwards thermodynamic theory, although of 
course it is far away of any physically feasible system. Nevertheless, it is interesting to 
qualitatively compare its predictions with the experimental results and with the behavior 
obtained by means of numerical simulation (Monte Carlo and molecular dynamics) of much 
more sophisticated and realistic models. 

The paper is organized as follows. In the next section the model to be used in the 
following is formulated as a Markov process, and some of its properties derived in ref. 17 1 
are shortly summarized. In particular, the steady distribution characterizing the long time 
limit of compaction is reminded and the Edwards compactivity identified. From the steady 
distribution, explicit expressions for some average properties are easily derived. Special 
attention is given to the relationship between the second moment of the volume fluctuations 
and the average volume, since it has been used to compute the compactivity in experiments. 
The detailed shape of the volume fluctuations is analyzed in Sec. IIIIl where it is predicted to 
be Gaussian in the limit of a large system, i.e. with a large enough number of particles. This 
is verified by means of Monte Carlo simulations of the model equations. The above model 
predictions are compared with experimental findings and simulation results for more realistic 
models in Sec. IIVI Special emphasis is put on the discussion of the boundary conditions, 
e.g. open versus closed systems, the local or global character of the measurements, and the 
concept of effective temperature and its relationship with the compactivity. In Sec. IVj the 
linear response of the model to a change in the tapping intensity is investigated. The volume 
response function is computed by Monte Carlo simulation, and a stretched exponential decay 
is found in the relevant part of the relaxation. Finally, in the last section, a short summary 
and some additional discussions are presented. 
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TABLE I: Elementary transitions and their transition rates for the lattice model. Particles and 
holes are represented by circles and crosses, respectively. 



Transition Initial state Final state Rate 

diffusion 00X0 0X00 p/2 

diffusion 0X00 00X0 p/2 

hole annihilation 0X0X0 0X00 p/2 

hole annihilation 0X0X0 00X0 p/2 

hole creation 0X00 0X0X0 q 

hole creation 00X0 0X0X0 q 



II. THE MODEL AND THE STEADY DISTRIBUTION 

Consider a one-dimensional lattice having + 1 particles and a variable number of sites. 
The dynamics of the system is modeled trying to capture, in a very siniplified way, the 
relevant features of compaction processes in tapped granular systems 0, y, Q- Then, a 
Markov process is introduced with the elementary transitions given in Table [H where the 



circles represent particles and the crosses empty sites or holes, and p and q are arbitrary 
positive real parameters . Time is measured in some arbitrary units proportional to the 
number of taps. Only those particles involved and/or conditioning each of the transitions 



are represented. It is worth to point out that the equa^ 



ity of the rates for diffusion and hole 



annihilation follows from the rational of the model jl7 |. 

The transitions in Table [T] define the modeled effective dynamics connecting metastable 
configurations in a system being tapped. These configurations are characterized by the 
property that all the holes are isolated, i.e. surrounded by two particles and, in order to 
enumerate them, a set of variables n = {ni,n2, . . . ,n^} is introduced. By definition, the 
variable rii takes the value = 1 if there is a hole to its right, separating it from particle 
i + 1. On the other hand, rij = if the site to the right of particle i is occupied by particle 
i + 1. It is assumed that there can not be a hole either to the left of particle 1 or to the 
right of particle + 1. Note that this feature is preserved by the dynamical events allowed 
in the model. Moreover, the dynamics does not change the relative order of the particles 
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in the lattice and, therefore, the identity of the particles plays no role in describing the 
configurations of the system. More details, as well as a discussion of the physical motivation 
of the model, are given in ref . l3] • 

Introduce RiU as denoting the configuration obtained from n by modifying the value of 
rii, keeping the same all the other variables, i.e.. 



RiH = {rii, nj_i, 1 - Tii, rii+i, um}- (1) 

Therefore, RiU differs from n only in the status, empty or occupied, of the site to the right 
of particle i. The probability transition rates W{n'\n) from configuration n to configuration 
n' defining the Markov process and represented in Table [H can be written in an analytical 
form as: 



1. Diffusion 



P 

W{RiRi+in\n) = - [(1 - 72^)72^+1 + ^^(1 - nj+i)] , (2) 



2. Hole annihilation 

P 

W{Rin\n) = - {ui^iUi + UiUi+i) , (3) 

3. Hole creation 

W{Rin\n) = q{l-ni){ni^i+ni+i). (4) 
This stochastic process can be proven to be irreducible. In addition, its (unique) steady 



probability distribution p*-'^-* (n) verifies the detailed balance condition and is given by 17|] 



?>'"(") = (S) 



where 7 = p/2g and 



N 



NH{n) = Y,ni (6) 



i=l 



is the total number of holes of the lattice configuration n. Finally, 2^(7) is the 'partition 
function' following from the normalization condition, 

N 

Z{l) = 5^7-^" = E ^(^h)!-''", (7) 

n Nh=1 
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with Q{Nh) being the number of metastable configurations having exactly Nh holes. A 
standard combinatorial reasoning yields 



From Eqs. and (IHl) it follows that 



Z(7) =[^ + -) - 1' (9) 
that for large N implies 

lnZ(7) ^ A^ln 1^1 + . (10) 

To be more precise, the required limit leading from Eq. (jH]) to Eq. fllOl) is ^ 1 and N ^ 'j, 
since if 7 is of the same order as iV 1, the correct limiting form for the partition function 
is 

Z(7) ^ e^/^ - 1. (11) 

In the following, attention will be restricted to systems having a number of particles such 
that Eq. flTUl) applies. 

Once lnZ(7) is known in a analytical form, it is a simple task to compute average values 
characterizing the macroscopic state of the system. Also, the probability distribution for 
the number of holes in the steady state, P^'^\Nh), is easily derived, 

Pi-\N„) ^ ^p")(n)^..,s„. = "'^l';:"" . (12) 

In the second equahty above, Sn,m is the Kronecker delta. The associated moment generating 
I I 

function |20j] is 

N ^ N 



Nh=1 ^^^^ Nh=1 



z{i) ■ 

Thus the cumulant generating function [20| of the variable Nh is 



(13) 



\nGN,{k)^N\n^-^^, (14) 
7 + 1 

where use has been made of Eq. ( fTOl) . From this expression, the steady average of the number 
of holes {Nh) and its variance a'j^ are directly computed. 



{NhY'^ = - 



9 In GN^{k) 
dk 

7 



a 



Nh 



(«)2 ^ _ 



(16) 



.=0 (7 + 1)^' 

From Eq. f|T5l) . it is seen that the steady state gets more compact as 7 increases. On 



the other hand, it is found that the relaxation of the system becomes slower 17|. This is 
consistent with the physical interpretation of the model, where is understood to play a 
role analogous to the tapping intensity in granular experiments j2|, ^ . Now, in order to get a 
canonical distribution of the form proposed by Edwards and coworkers 31], the 'compactivity' 
X is defined by 

X = 



1 



In 7 

so the steady distribution in Eq. ([5]) can be rewritten as 

-Nh/X 



ZiX) 



Z(X) = (l + e-i/^) 



N 



(17) 



(18) 



In terms of the compactivity, Eqs. f|T5l) and f[T6l) take the form 



(Nh) 



cr 



eiA + r (eV^ + l) 

respectively. Combination of the above two equations gives 



2 ' 



(19) 



,1/2X 



(20) 



showing the typical N'^/'^ dependence of the relative fluctuations of an extensive quantity 
on the number of particles of the system. 

Since p and q can take arbitrary positive real values, it is clear from the definition of 7 
given below Eq. ([S]), that its range of variation is < 7 < 00. Consequently, X can take 
any real value, i.e. —00 < X < +00. More precisely, when 7 is increased monotonically 
from to 00, the compactivity X decreases from 0^ to —00 (for 7 = 1~), then jumping to 
+00 (7 = 1"*") and decreasing again towards 0"^. Then, there is a parameter region, namely 
p < 2q, or < 7 < 1, in which negative values of the compactivity occur. With reference 
to the dynamics of the model, this happens when the transition rate for the hole creation 
processes is larger than that for hole annihilation processes. 

Equation ([T5]) lead to the relation 

d 



d{l/X) 



(Nh) 



(21) 
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and considering a series of experiments with the same number of particles in which the 
average number of holes is changed between {Nh)^i^ and {Nh)2\ 

(^h)^^' 1 1 1 

Here Xi and X2 are the initial and final values of the compactivity, respectively. The 
generalization of this equation for three-dimensional systems in the context of the Edwards 
theory was first noted by Nowak et al. [4], and used in ref. 6j] to obtain the compactivity as 
a function of the volume fraction from experimental data of a system submitted to pulses. 
More will be said about this in Sec. IIV[ 



III. VOLUME FLUCTUATIONS 

The length L of the system, measured in units of the distance between consecutive sites, 

is 

L = N + Nh. (23) 
In the steady state, the average value of L is 

(L)(^) = N+ {Nh^^^ = ^^t±^_ (24) 

7 + 1 

It is also convenient to introduce the length per particle I = L/N, so (l)^'^ = (7 + 2)7(7+1), 
or 

This can be transformed into a relationship between X and {l)^^^ by means of Eq. ([1 



1 2 - 

It is seen that 1/X is a monotonically decreasing function of vanishing for {l)^^^ = 

3/2 that, of course, corresponds to 7 = 1. Therefore, it is concluded that X > for 
(/)*^*) < 3/2, while X < for (Z)'-'*^ > 3/2. Negative compactivities correspond to the loosest 



configurations and have also been obtained in other models of compaction [2ll]. The inverse 
of the compactivity does not vanish in the limit of the loosest possible configuration, but it 
tends to —00. Obviously, the existence of negative values of the compactivity does not violate 
any general physical principle and does not imply any instability of the steady state. They 



are a consequence of the fact that the number of configurations fl{NH) is not a monotonic 
increasing function of the length of the system (or Nh), but it exhibits a maximum, going 
afterwards to unity for the maximum value of the number of holes, = N. Equivalently, 
the possibility of steady states with X < is directly related with the length L of the 
system having a finite upper bound. From a formal point of view, the situation is similar to 
the existence of negative temperatures in systems of nuclear spins in certain paramagnetic 
crystals placed in an external magnetic field. 

From Eq. (fT6ll . the variance of the total length can be written as 

^i = < = iV((0^^^-l)(2-(0^^^)- (27) 

In order to get a more detailed information of the length fluctuations, introduce the stochas- 
tic variable Y by 

The moment generating function of Y is 

Gy{k) = (e'^^)(«) = e {e^^H)(') 

= e ^-n (29) 

where Gjq^ was defined in Eq. (fT3l) . Then, using Eqs. (fT4|) - (fT6|l . 

InGy(fc) = -2fcAri/27-V2 + jvinjl + 7-iexp [iN-^/^-i-^''^{l + -i)k]] - Nln^^^, (30) 

that in the limit ^ 1 reduces to 

lnGy(A;)~-^. (31) 

This proves that, for a large enough number of particles, the stochastic variable Y obeys 
a Gaussian distribution with vanishing average value and unit variance. Formally, this is 
equivalent to saying that the probability density uji{l) for the variable / is also Gaussian with 
average (Z)'-'*^ and standard deviation ai = ol/N. 

A numerical check of the Gaussian character of Y for large values of is presented in 
Fig. [H where Monte Carlo simulation results for the steady probability distribution P^^^ {Y) 
are shown for several values of 7. The number of particles used is A^ = 10^. The solid line 
is the Gaussian distribution of zero average and unit variance. An excellent agreement is 
observed, even for rather small values of P^^^iY). 
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FIG. 1: Steady probability distribution of the dimensionless variable Y defined in Eq. (j28p . 

The symbols are Monte Carlo simulation results obtained with different values of 7, as indicated 
in the insert. The solid line is the Gaussian theoretical prediction. 



In order to compare the results derived for the model in this Section with some experi- 
mental findings, it is convenient to consider also the volume fraction cj) defined as = 
As a consequence of cx; being proportional to A^~^/^, the probability density for 0, is 
also Gaussian, again in the limit of large A^, with average and variance given by 



and 



[((/)(-) -1) (2- (/)(-))] 



respectively. Thus Eq. f l26l) is, for large A^, equivalent to 



X 1 - 



(32) 



(33) 



(34) 



IV. COMPARISON OF THE MODEL PREDICTIONS WITH SOME EXPERI- 
MENTAL RESULTS 



Of course, given the simplicity of the model discussed here and its one-dimensional char- 
acter, it can not be expected to quantitatively reproduce the experimental and simulation 
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results. Nevertheless, it is interesting and illuminating to qualitatively analyze those results 
in the light of the model predictions. 

The first experimental investigation of the volume fluctuations in the steady state of a 
tapped granular medium seems to be the one reported in ref. j4]. More exactly, the property 
considered in that work is the density in a given region, located at a certain height of the 
system. Then, the volume is kept constant and the measured density changes are due to 
variations in the number of particles within the volume considered. The results indicated 
that the probability distribution of the density fluctuations exhibited in the majority of 
the cases a Gaussian shape, as found in the model, although some significant non-Gaussian 
deviations were also found, particularly near the bottom of the vibrated column. Moreover, 
the dependence on the number of particles of the relative variance was also the same as in 
the model, i.e N~^^'^ as shown in Eq. fl20|) . 

Variations in the inverse of the compactivity were also measured, by using the relationship 
between the second moment of the volume fluctuations ay and the compactivity following 
from Edwards' theory, i.e. the generalization of Eq. (1221) . To do so, the authors considered 
the specific volume v defined as the inverse of the packing fraction. This corresponds to the 
quantity / in our model. The reported results are 

al = a + h{v)\ (35) 

^ oc ln(a + + constant, (36) 

j\ 

where a and h are two dimensionless parameters. Taking into account that the available 
range of experimental data is rather limited and located in the proximity of the close packing 
limit, it can be fairly concluded that the above expressions are in qualitative agreement with 
the predictions of the model, namely the expressions for o"/ and X~^, Eqs. fl34p and fl26p 
respectively. Nevertheless, to put this comparison in a fair proper context, the following two 
comments seem appropriated: 

1. In ref. Q], the Edwards theory is considered in a local form, i.e. it is assumed that 
the state of each subregion of the system in which measurements are carried out 
is described by a canonical distribution with a different value of the compactivity. 
Although this might be a sensible idea, it is clearly an extension of the theory as 
originally formulated by Edwards. 
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2. Moreover, what is measured in the reported experiments are the fluctuations in the 
number of particles at constant volume. This is different from the fluctuations in the 
volume at constant number of particles, although both can be associated to density 
fluctuations. In fact, it can be shown that the quantity actually obtained from the 
second moment of the number of particles fluctuations by means of a relationship 
similar to Eq. (!22|) is not the compactivity, but some kind of 'fugacity' [ITj. Although 
both quantities are simply related and perhaps their values are very close in real 
granular systems, they are conceptually quite different. 

Another experimental study of the steady volume fluctuations in a granular medium has 
been carried out by Schroter et al. Q]. This work might look as quite similar to the one 
discussed above, but they differ in important aspects, both in the methodology and in the 
results. In ref . [Gj , the fluctuations of the total volume of the system were measured, although 
the results are expressed in terms of the volume fraction (p. A Gaussian distribution was 
found in all the reported cases. Also, the ratio between the standard deviation and the 
average turned out to be proportional to the inverse of the square root of the number 
of particles. Both results are the same as derived above for the model considered here. 
Nevertheless, the variance of the volume fraction fluctuations did not show the analogous 
of the simple linear behavior reported in ref. that, as it has been already mentioned, 
is consistent with the linear behavior predicted by Eq. (!34|) near the close packing limit. 
Instead, it presents a well deflned minimum. A flrst, naive explanation is provided by the 
fact that a wider density range is analyzed in than in 4]. Then, it could be expected that 
the linear law in the latter corresponds to an approximate description of a small window of 
the density interval considered in the former, namely to one with the largest volume fraction 
values. Even assuming this explanation, the origin of the minimum would remain as an open 
question. 

The result obtained with the model, Eq. (!34|) . when considered over the density range 
1/2 < < 1, also shows a non-monotonic behavior, but presenting a maximum instead 

of a minimum exhibited by the experimental data in ref. 6j. Then, it seems plausible to 
conclude that the existence of the minimum, if conflrmed, must be due to some effects that 
are not captured by the simple model considered here. Given that it appears for high values 
of the density, it is tempting to speculate whether the increase of the volume fluctuations in 
that region, to the right of the minimum, is an indication of the presence of some 'singular' 
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behavior associated, for instance, with the existence of a maximum random packing, a 
phenomenon that is not present in the modeL 

The Edwards compactivity X was also determined by Schroter et al. by means of the 
granular version of the equilibrium fluctuation-dissipation theorem. To get explicit values 
for X, these authors assume that the inverse of the compactivity vanishes in the random 
loose packing limit. Since X is a monotonic decreasing function of the volume fraction, 
this boundary condition implies that it is always positive. As it has already been indicated, 
there is nothing in the granular statistical mechanics theory requiring this to be true. If 
attention is restricted to the parameter region in which the model has positive compactivity, 
i.e. (0)*^'*-' > 2/3, the shape of the curve X{{(j))^^^)^ Eq. ([31]), is similar to the one reported 
in Fig. 4(c) of ref. Q]. 

Although the recent results presented in ref. g] have not been obtained experimentally, 
but by means of molecular dynamics and Monte Carlo simulations, they deserve special 
attention because of their relevance and the interest of the analysis carried out. In particular, 
interest is focused on a parameter defined as the 'configurational granular temperature', 
instead of the compactivity, although some close relationship between both parameters is 
claimed to hold. A model system of grains with normal and tangential forces is considered. 
The instantaneous volume fraction is measured in the steady state and analyzed as a function 
of its average value. To measure them, a region in the bulk of the system is considered. It 
is important to stress that this region is defined in such a way that both its volume and its 
number of particles fluctuate in time. Consistently with all the results discussed above and 
the model predictions, Gaussian fluctuations of the volume fraction are also found in the 
simulations. 

The analysis of the results in 8| differs, at least conceptually, of the thermodynamic 
theory as proposed initially by Edwards [3]. Instead of considering that the macroscopic 
parameter characterizing the state of the system is the volume, the energy is used like in 

is] 

normal molecular systems. More precisely, the probability pr of finding the system in the 



blocked metastable state r of energy Er is assumed to be given by [10|, ll8| 

oc e-^'-Z^-f, (37) 
where Tconf is referred to as the configurational granular temperature of the system. The 
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above expression is to be compared with Edwards' proposal, see also Eq. (fT8|) . 

pis) ^ g-K/x^ (3g) 

Vr being the volume of the system in the configuration r. Both expressions only become 
equivalent if, for the relevant configurations, i.e. those with a non-vanishing probability. 



the energy and the volume are proportional 



22j, with a configuration-independent constant. 



Consequently, the configurational temperature Tconf and the compactivity X also would be 
proportional in this case. 

In the paper, the authors plot (fig. 4) the standard deviation o"^ of the volume fraction 
fiuctuations as a function of the steady average volume fraction {(f))^^^. In the plotted interval, 
a linear behavior of a^{{(f))^^^) is clearly identified. Then, by usin g th e relationship between 



afj) and the compactivity X, an analytical expression of the form 23| 

^ = c + /((0)W) (39) 

is obtained. Here, the explicit form of the function f{{(f))^'^^) is known, while c is an unknown 
constant. If Tconf = dX, with d a constant, Eq. fl39|) yields 

fm^'^) = 7f^-c. (40) 

conf 

In ref. [8], a plot of the volume fraction as a function of the configurational temperature, 
obtained by means of Monte Carlo simulation, is also provided (fig. 5). Consistently, this 
plot covers the same density range as the one giving the fiuctuations. Therefore, it is posible 
to check whether the linearity between f{{(f))^^^) and T^~nf predicted by Eq. pUj) is verified. 
When this is done, a systematic but quite small deviation from linearity is observed. Then, 
taking into account possible numerical errors as well as the fact that the measurements are 
carried out in such a way that both the volume considered and the number of particles 
inside it fluctuate, no definite conclusion can be reached on the equivalence of the canonical 



distribution in energy, Eq. (1371) . and volume, Eq. (138!) [231 ]. This is an important point. 



both from fundamental and applied perspectives and deserves further attention. 



V. LINEAR RESPONSE. VOLUME AUTOCORRELATION FUNCTION 



In this Section, the linear response of the model to a perturbation of the parameter 7, 
introduced below Eq. ([5]) will be studied. This parameter measures the ratio of the transition 
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rates of the processes associted with hole annihilation and creation, respectively. Because 
of the definition in Eq. f|T7|) . a change in 7 is equivalent to modifying the compactivity X. 
In the context of real experiments, it might correspond to varying the vibration intensity, 
defined as the ratio between the peak acceleration of the shakes and the gravity acceleration 

0. 

The typical relaxation experiment to be considered is as follows. The system is initially 
prepared in a uniform density configuration corresponding to a low density state, namely 
with Nh = N. Then, it is allowed to evolve with a value 7 + A7 of the rate ratio, until it 
reaches a steady state with a probability distribution p^^\n; 7 + A7). At a given moment, 
taken as the time origin t = 0, the rate ratio is instantaneously changed to the value 7. It is 
assumed that IA7I <^ 7. Then, the relaxation of the system to a new steady state j9*^*^(n; 7) 
is followed. 

The probability distribution function of the system at t = is 

p(n, t = 0) = p(^) (n; 7 + A7) ^ p(^) (n; 7) + ^&^A-f. (41) 

c'7 



(42) 



Using Eq. (]T8|) . it is obtained 



apW(n;7) _ gpW(n;7) rfX _ iVj, - (jV^; 7)^^^ .dX 
d-i dX rf7 X2 ^ ^"■'^^^7 



The average value of a function F{n) of the system configuration at time t > is, by 
definition, 

(F,t) = 5^F(n)p(n,t); (43) 

n 

where the probability distribution p{n, t) is obtained from p{n, t = 0) through 

p{n,t) = ^p{n,t\n',0)p{n',t = 0). (44) 

n' 

Here, p{n, t\n', 0) is the conditional probability that the configuration of the system is n at 
time t given it was n' at t = 0. By construction, p^^\n]'y) is stationary for the dynamics 
occurring for t > 0, i.e. with rate ratio 7, 

J]p(n,t|n',0)pW(n';7) =pW(n;7). (45) 

n' 

Taking this into account and using Eqs. fHTl) and fH2l) . Eq. (l43l) becomes 

1 dX 

{F, t) = {F; 7)(^) + (AF(t) AL(0); 7)^^) A7, (46) 
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where (AF(t)AL(O); 7)*^^^ is the steady time correlation function of F{n) and the volume L 
defined as 

(AF(t)AL(0);7)(^) = ^ ^ AF(n)p(n, t|n', 0)AL(n')p('H^'; 7), (47) 

n n' 

with 

AF(n) = F(n) - (F;7)(^\ (48) 

AL(n) = L(n) - (L;7)(^). (49) 

The normalized relaxation function for the property F, $ir(t) is 

^ _ {F;t)-{F;^y^^ _ (AF(t)AL(O); 7)^-) 
"^^^ = (F,0)-(F;7)^ " (AFAL;7)W ' ^ ^ 



For the special choice F = L, 



^ <^"W^^'»)^^>"'. (51) 



This expression is an exact consequence of modelling the dynamics of the system by means 
of a Markov process and of the canonical form of the steady probability distribution. The 
Monte Carlo simulation results for the relaxation function to be presented in the following, 
have been obtained by using Eq. fl^ . Figures [31 and S] show $L(i) for three values 
of the rate ratio, namely 7 = 5, 15, and 50, respectively. Other values have also been 
investigated, and the results are consistent with the comments carried out below. In the 
figures the quantity ln(— ln$j^) has been plotted as a function of Int. It turns out that 
this is a convenient representation to identify the shape of the decay of $L(t). Note that 
in this representation, a straight line corresponds to a stretched exponential or Kohlrausch- 
Willians- Watts (KWW) function, 

^Kwwit) = exp 



(52) 



where r and (3 are the two parameters identifying the function. 

In the three figures, two different regions are easily identified. For 'short' times, there is 
a linear region with slope /5 ~ 1. i.e. an exponential decay. Afterwards, there is a region of 
'intermediate' times in which another linear behavior shows up, but now with (3 < 1. The 
initial region extends over a time interval, Int < 2, that is roughly independent of 7. In 
this time interval, the response function $i decays a small amount. Moreover, this initial 
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FIG. 2: Double logarithm plot of the normalized response function for the volume as a function 
of the logarithm of the time for 7 = 5. The latter is measured in dimensionless units proportional 
to the number of Taps. The dashed line is the Monte Carlo simulation result, while the two solid 
lines are linear fits as discussed in the text. 

decay becomes smaller and smaller as the value of 7 increases. On the other hand, in the 
second time window in which the response function shows a KWW behavior, most of the 
relaxation occurs. For 7 = 5, Fig. ([2]), the fit to a KWW function with (3 ~ 0.659 is quite 
good for the time interval 3 < Int < 6, that corresponds to —2 < ln(— ln$x) < 0.75, i.e. 
0.87 ^ $L ^ 0.12. For 7 = 15, Fig. [HI the value of (3 in the intermediate time window is 
(3 ~ 0.567, and the interval of the relaxation function described by it is 0.92 ^ $l ^ 0.12. 
Finally, for 7 = 50, Fig. H it is /5 ~ 0.534 for 0.95 > > 0.12. 

The kind of behavior of the response function has been often found in the context of 
structural glasses, both in experiments and in simple models, for the linear relaxation of the 



energy following a temperature perturbation 2% l25|, |26|| . Then, this is another indication 



of a possible conceptual connection between structural glasses and dense granular systems. 
Moreover, the simulation results for the model clearly indicate that the value of the exponent 
(3 characterizing most of the relaxation tends to 1/2 as the value of 7 increases. In simple 



26 



27| 



models of structural glasses, this value is associated to diffusive processes 

Although an explicit analytical solution of the model has not been derived yet, it is 
possible to understand the diffusive origin of the exponent [3 = 1/2 as well as to get some 
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FIG. 3: The same as Fig. [21 but for 7 = 15. 
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FIG. 4: The same as Fig. [21 but for 7 = 50. 

information about the characteristic relaxation time r. For 7^1 [X ^ 1), a typical 
metastable configuration of the system consists of isolated holes separated a distance of the 
order of 7 (see, for instance, Eq. (IT^ ). In this situation, the only relevant process to decrease 
the volume is diffusion: holes move through the system until two of them get together and 
one is destroyed. In this picture, the characteristic relaxation time will have the diffusion 
form r oc 7^ = e^^-^ ^ i.e. an Arrhenius type law with the molecular temperature replaced by 
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FIG. 5: Plot of the relaxation time r appearing in the KWW expression, Eq. (I52p . as a function 
of the rate ratio 7. Time is measured in units proportional to the number of taps. The symbols 
are from the Monte Carlo simulation data and the solid line a linear fit of them. 

X. 

In Fig. [5], the Monte Carlo simulation results for r, obtained by fitting the relevant 
intermediate time window to a KWW function as discussed above, are plotted as a function 
of 7. Note the logarithmic scale, so that the observed straight line represents a power law 
behavior. The best fit, given by the solid line, yields r 6.877^'°^ in good agreement with 
the result following from the diffusive picture. 

VI. CONCLUSION 

A simple model for compaction in granular media has been used to investigate some 
static and dynamic properties of the steady sate reached by the system in the long time 
limit. The simplicity of the model allows to obtain analytical expressions for some of its 
properties. On the other hand, it reproduces qualitatively well many of the properties found 
in real granular systems. In the context of static properties, this includes the Gaussian 
shape of the volume fluctuations, the dependence of its variance on the average volume and 
the number of particles, and the possibility of identifying the compactivity as defined by 
Edwards. 
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Moreover, the scenario obtained from the analysis of the model, provides an useful tool 
to interpret the experimental results and molecular djTiamic simulations of more realistic 
models. This is exemplified in the critical revision carried out in Sec. \IV\ where the need 
of a clear differentiation between density fluctuations at constant volume and at constant 
number of particles shows up. Also, attention must be paid to distinguish between local and 
global applications of the thermodynamic theory. 

Another issue that is worth commenting is the possible existence of negative values of 
the compactivity and/or the effective temperature. If the mechanical statistical description 
is assumed as the starting point of a thermodynamic description, there is no solid reason to 
deny by principle this possibility. Even more, one is tempted to conclude that it is highly 
probable, at least at a theoretical level. Negatives values of the compactivity will show 
up unavoidably if the number of metastable configurations is not a monotonic increasing 
function of the volume. The possible volume occupied by the tapped granular medium is 
limited by the loose and random jamming packings. In both limits, the number of possible 
metastable configurations can be expected to be quite reduced. The requirement of maxi- 
mum or minimum space between particles impose a rather severe restriction on the allowed 
arrangement of the particles. Therefore, it seems sensible to expect that the number of 
configurations for intermediate values of the volume be larger than at the extremes, present- 
ing a maximum at some values. The same qualitative reasoning applies to the effective or 
configurational temperature changing the volume by the energy. 

The question then is why negative compactivities (effective temperatures) have been not 
identified in the experiments yet. There are two main reasons for it. The first one is that 
only differences of compactivities have been measured up to now, and it has often been 
arbitrarily assumed that the compactivity diverges in the random loose packing limit. The 
second reason is that negative values of the compactivity, if they exist, would correspond to 
the lowest density theoretically accesible region, and it can be very hard to reach in practice. 
In any case, this is a point that definitely deserves much more study. 

In the last part of the paper, the linear relaxation of the model has been considered. A 
slow relaxation, accurately represented by a KWW function with exponent (3 tending to 
1/2 as the parameter corresponding to the vibration intensity is decreased, has been found. 
Moreover, the expression of the characteristic relaxation time r has Arrhenius form (with 
the consistent substitution of the compactivity by the temperature). 
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The behavior of the hnear response function for the model is fully similar to what is 
observed in the relaxation of many structural glasses and also in models of glassy relaxation. 
It must be noted that the model discussed here was built to mimic the compaction experi- 
ment and that, consistently, the relaxation of the density starting from a loose configuration 
exhibits the characteristic inverse logarithmic law [l^ found in experiments {3, 4|. It seems 
interesting to experimentally investigate whether linear response in tapped granular systems 
is also slow and described by a KWW function. If that were the case, it would be another 
indication of similarity between two apparently different kind of systems, granular media 



and structural glasses, reinforcing the link already found in other phenomena 281], including 



the existence of hysteresis cycles {4, |29| and memory effects 
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3l|. 



Acknowledgments 



This research was supported by the Ministerio de Educacion y Ciencia (Spain) through 
Grant No. FIS2005-01398 (partially financed by FEDER funds). 



[1] G.C. Parker and A. Mehta, Phys. Rev. A 45, 3435 (1993). 

[2] J.B. Knight, C.G. Fandrich, C.N. Lau, H.M. Jaeger, and S.R. Nagel, Phys. Rev. E 51, 3957 
(1995). 

[3] S. F. Edwards and R. B. S. Oakeshott, Physica A 157, 1080 (1989); S. F. Edwards and A. 
Metha, Journal de Physique 50, 2489 (1989); S. F. Edwards and C. C. Mounfield, Physica A 
210, 279 (1994). 

[4] E.R. Nowak, J.B. Knight, E. Ben-Nairn, H.M Jaeger, and S.R. Nagel, Phys. Rev. E 57, 1971 
(1998). 

[5] P. Philippe and D. Bideau, Europhys. Lett. 60, 677 (2002). 

[6] M. Schroter, D.J. Goldman, and H.L. Swinney, Phys. Rev. E 71, 030301(R) (2005). 
[7] H. A. Makse and J. Kurchan, Nature 415, 614 (2002). 

[8] M.P. Ciamarra, A. Coniglio, and M. Nicodemi, Phys. Rev. Lett. 97, 158001 (2006). 
[9] A. Barrat, J. Kurchan, V. Loreto, and M. Sellitto, Phys. Rev. Lett. 85, 5034 (2000); A. 
Barrat, J. Kurchan, V. Loreto, and M. Sellitto, Phys. Rev. E 63, 051301 (2001); V. Cohzza, 



22 



A. Barrat, and V. Loreto, Phys. Rev. E 65, 050301 (2002). 
[10] A. Fierro, M. Nicodemi and A. Coniglio, Europhys. Lett. 59, 642 (2002). 
[11] D. S. Dean and A. Lefevre, Phys. Rev. Lett. 90, 198301 (2003). 

[12] J. J. Brey, A. Prados, and B. Sanchez-Rey, Phys. Rev. E 60, 5685 (1999); Physica A 275, 310 
(2000). 

[13] A. Lefevre and D. S. Dean, J. Phys. A 34, L213 (2001). 

[14] J. Berg, S. Franz, and M. Selhto, Eur. Phys. J. B 26, 349 (2002). 

[15] G. De Smedt, C. Godreche, and J. M. Luck, Eur. Phys. J. B 27, 363 (2002). 

[16] A. Prados and J. J. Brey, Phys. Rev. E 66, 041308 (2002). 

[17] J.J. Brey and A. Prados, Phys. Rev. E 68, 051302 (2003). 

[18] A. Conigho and M. Nicodemi, Physica A 296, 451 (2001). 

[19] E.R. Nowak, J.B. Knight, M. Povinelh, H.M. Jaeger, and S.R. Nagel, Powder Technol. 94, 79 
(1997). 

[20] N.G. van Kampen Stochastic Processes in Physics and Chemistry (North-Holland, Amster- 
dam, 1992). 

[21] R. Monasson and O. Pouhquen, Physica A 236, 395 (1997). 

[22] EPAPS Document No. E-PRLTAO-97-045642. For information on EPAPS, see 

|http:/ /www.aip.org/pubservs/epaps.html 
[23] J.J. Brey and A. Prados, unpublished. 

[24] See, for instance, G.W. Schrerer, Relaxation in Glass and Composites (Wiley, New York, 
1986). 

[25] G.H. Fredrickson, Annu. Rev. Phys. Chem. 39, 149 (1988). 
[26] J.J. Brey and A. Prados, Phys. Rev. E 53, 458 (1996). 

[27] J.J. Brey, A. Prados, and M.J. Ruiz-Montero, J. Non-Cryst. Solids 172-174, 371 (1994). 
[28] D.I. Goldman and H.L. Swinney, Phys. Rev. Lett. 96, 145702 (2006). 
[29] A. Prados, J.J. Brey, and B. Sanchez-Rey, Physica A 284, 277 (2000). 

[30] C. Josserand, A. Tkachenko, D.M. Mueth, and H.M. Jaeger, Phys. Rev. Lett. 85, 3632 (2000). 
[31] J.J. Brey and A. Prados, Europhys. Lett. 57, 171 (2002). 



23 



